/*------------------<--- Start of Description -->--------------------\
| Power for comparison of two median survival curves;                |
| Provides approximate power for logrank test.                       |
| Assumptions are uniform accrual and exponential survival.          |
| Arguments:                                                         |
|   - Required:                                                      |
|     T0 = length of accrual period(uniform accrual assumed)         |
|     T  = total length of study(=T0+follow-up after last patient)   |
|     n1,n2 = sample sizes for groups 1 and 2 respectively           |
|     y1 = true group 1 median survival                              |
|     min_y2 = smallest possible true group 2 median survival        |
|   - Optional:                                                      |
|     max_y2 = largest possible true median survival for group 2     |
|     inc_y2 = increment value for range of y2                       |
|     alpha = type 1 error, e.g. .01 or .05, default=.05             |
|     sides = 1 or 2 for 1 or 2 sided test, default=2                |
|     plot  = 'P' for line printer plot of group 2 sample size vs y2 |
|             'G' for SAS/GRAPH plot of group 2 sample size vs y2    |
|     unit   = units for time, e.g. years, months, days, hours       |
| Output: Power for true group 1 mean of y1 vs true group 2 means    |
|         ranging from min_y2 to max_y2                              |
|--------------------------------------------------------------------|
|----------------<--- Start of Example and Usage -->-----------------|
| Example: %sv2md_pr(sides=1,n1=188,n2=188,y1=3.466,min_y2=0.31,     |
|                    max_y2=2.31, inc_y2=.2, t0=3,t=5,unit=yrs,      |
|                    plot=p);                                        |
| Usage: %sv2md_pr (ALPHA=.05,SIDES=2,N1=.,N2=.,Y1=.,MIN_Y2=.,       |
|      MAX_Y2=.,INC_Y2=.,T0=.,T=.,PLOT= , UNIT=);                    |
| Reference: Bergstralh, EJ.  SAS macros for sample size and power   |
|            calculations.  Proceedings of the 9th annual SAS Users  |
|            Group International Conference.                         |
|            Equation #18.                                           |
\-------------------<--- End of Example and Usage -->---------------*/
%MACRO sv2md_pr (ALPHA=.05,SIDES=2,N1=.,N2=.,Y1=.,MIN_Y2=.,
       MAX_Y2=.,INC_Y2=.,T0=.,T=.,PLOT= , UNIT=);
/*--------------------------------------------\
| Author:  Michael Riggs and Eric Bergstralh; |
| Purpose: Power for comparison of two median |
|          survival curves;                   |
\--------------------------------------------*/
 OPTIONS MISSING=' ' NOCENTER;
   %LET PLOT=%UPCASE(&PLOT);
 DATA T1;
      ALPHA=&ALPHA;
      SIDES=&SIDES;
      Y1=&Y1;
      MIN_Y2=&MIN_Y2;
      MAX_Y2=&MAX_Y2;
      INC_Y2=&INC_Y2;
      N1=&N1;
      N2=&N2;
      T0=&T0;
      T=&T;
    ZALPHA=(PROBIT(1-ALPHA))*(SIDES=1) + (PROBIT(1-ALPHA/2))*(SIDES=2);
    IF MAX_Y2=. THEN DO;
      MAX_Y2=MIN_Y2+1; INC_Y2=MIN_Y2+2;  *NEED 1 EXEC OF DO;
    END;
    TY1=Y1;
    TT=T;
    TT0=T0;
    TN1=N1;
    TN2=N2;
    CALL SYMPUT ('TOTFOL', PUT((T-T0),4.1));
    DO Y2=MIN_Y2 TO MAX_Y2 BY INC_Y2;
       LAMBDA1=LOG(2)/TY1;
       LAMBDA2=LOG(2)/Y2;
       LAMBD_BR=(TN1*LAMBDA1+TN2*LAMBDA2)/(TN1+TN2);
       IF TT0 GT 0 THEN DO;
         PHILBDA1=LAMBDA1**2   *
                  1/ (1-(EXP(-LAMBDA1*(TT-TT0))-EXP(-LAMBDA1*TT)) /
                               (LAMBDA1*TT0) );
         PHILBDA2=LAMBDA2**2   *
                  1/ (1-(EXP(-LAMBDA2*(TT-TT0))-EXP(-LAMBDA2*TT)) /
                               (LAMBDA2*TT0) );
         PHILBDBR=LAMBD_BR**2  *
                  1/ (1-(EXP(-LAMBD_BR*(TT-TT0))-EXP(-LAMBD_BR*TT)) /
                               (LAMBD_BR*TT0) );
       END;
       IF TT0=0 THEN DO; *ALL PTS ENTER STUDY AT SAME TIME;
          PHILBDA1=LAMBDA1**2/(1-EXP(-LAMBDA1*TT));
          PHILBDA2=LAMBDA2**2/(1-EXP(-LAMBDA2*TT));
          PHILBDBR=LAMBD_BR**2/(1-EXP(-LAMBD_BR*TT));
       END;
       ZBETA=(ABS(LAMBDA1-LAMBDA2)-
               ZALPHA*SQRT(PHILBDBR*(1/TN1+1/TN2)))/
             SQRT(PHILBDA1/TN1+PHILBDA2/TN2) ;
       POWER=PROBNORM(ZBETA);
       OUTPUT;
     END;

LABEL Y1="Group 1*Median (&UNIT)*Survival"
       Y2="Group 2*Median (&UNIT)*Survival"
       Power="Power";
 PROC PRINT SPLIT='*';
  ID  Y1; var Y2 POWER;
 TITLE3
   "POWER ESTIMATES FOR COMPARISON OF TWO MEDIAN SURVIVALS" ;
 TITLE4"Uniform accrual for &t0 &unit(T0) is assumed with analysis at &t &unit
 (T)";

 title5
 "Alpha=&alpha, Sides=&sides, N1=&n1, N2=&n2, True Group 1 Median Survival=&y1
  &unit";

 %IF &MAX_Y2 NE . %THEN %DO;
 %IF &PLOT= P  %THEN %DO;
       PROC PLOT NOLEGEND;
           PLOT POWER*Y2/ VAXIS=0 TO 1.0 BY .10
                          HAXIS=&MIN_Y2 TO &MAX_Y2 BY &INC_Y2;
      Label     Y2="Group2 median survival(&unit)";

    %END;
    %ELSE %IF &PLOT= G  %THEN %DO;
          PROC GPLOT ;
      PLOT POWER*Y2/ vaxis=0 to 1 by .2;
      SYMBOL1  F=SPECIAL V=J H=1 I=j       L=1;
      SYMBOL2  F=SPECIAL V=M H=1   I=j       L=2;
      Label     Y2="Group2 median survival(&unit)";
    %END;
    RUN;
%END;
%MEND sv2md_pr;